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Work during the period 

During the past six months effort was directed to the 
following areas: 

1. Procurement and installation of the line printer, 

% 

A-D converter boards, clock and extra memory for the MINC 11-23 
minicomputer, and the logic analyzer, 

2. Interconnection of the computer inputs to magnetic 
balance outputs to provide computer position control and data 
acquisition. This necessitated construction of four junction 
boxes and two 27 pair cables 70 feet long to provide connec- 
tions from balance to computer and computer to balance. 

3. Completion of the investigation of the use of parameter 

identification of a means of determining dyncimic characteristics. 
This study indicated that the technique was successful at pro- 
viding preliminary measurements of for the ogive 

a *"q 

cylinder model. Sufficient measurements were not made to pro- 
vide a comparison with other methods; i.e., forced oscillation 
motion. The measurements made indicated that results could be 

obtained with random noise inputs of .1, .2 and .4 degrees in 

* 

pitch. This is reported in Gautham Ramohalli's M.S. thesis , 
an abbreviated version of which constitutes the bulk of this 
report. 

4. Repair of Thyratron and motor generator power supplies 
for the pitch and yaw degrees of freedom. Vacuum tubes and de- 
fective wiring in the ThyratrvOn power supply used in pitch were 

Stability Derivative Estimation using Digital Signal Processing 
Methods with the Magnetic Suspension Wind Tunnel, S.M. Thesis, 
Dept of Aeronautics and Astronautics, MIT, May, 1981. 
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replaced to permit use of the supply during extensive repair 
of the motor generators. Time since overhaul on these supplies 
was about 1100 hours, which is 100 hours longer than the expected 
life of the generator bearings. 

Because of the increased effort required for computer 
connection and system repairs, grant funds have been expended 
as of this reporting period. 

Personnel 

During the past six months the following professional 
personnel have contributed to the effort under this program i 

Professor Eugene E. Covert 
Dr. Charles W. Haldeman 
Mr. Gauthcbn Ramohalli 

Choice of a method for 
handling dynamic system data 

The first method considered was the EEIC method of Merhav 
and Gabay (1) . This has the attractive feature of being able 
to estimate not only poles but zeroes of a continuous system 
directly. Further, the method drives the pole zero excess co- 
efficients (if a higher order model is guessed) automatically 
to zero and any excess pole zero pairs are driven to each other. 

This method was tried on a simple first order model on a 
simulation on the VAX 11/780. The method showed good con- 
vergence to reasonably accurate values. However, it needed 
about 50 seconds of data before the values settled down. As 
the order of the mod'll is increased, the complexity of the 
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process increases. For each order increase an additional 

V 

filter is required of the form s/(s+a) . Further, the pole 
zero pairs need to be monitored and factorization used to 
cancel any poles and zeroes that converge to each other (but 
do not converge to a steady valuel). This method was not used. 

If the memory and computational requirements are to be 
minimized, analog filters are indicated. However, since 
evaluating digital data processing for the application is an 
important goal, the data was processed using the M/NC. 

Since the MINC was used to acquire many data points in 
a short time, it was located near the wind tunnel. To permit 
its use for data reduction as well as acquisition, a reduction 
method was selected that could be easily implemented on the 
MINC rather than having to transfer the data over to the 
larger central VAX system. This led to choice of the second 
method . 

The second method considered was that of the Fourier 
spectrum analysis using Fast Fourier Transforms algorithms. 
Using the Institute of Electrical and Electronics Engineers 
digital signal processing program (2) , FFT subroutines were 
implemented on the MINC 11/23 described below. However, 
spectrum analysis methods have the shortcoming that their 
resolution in fivequency is limited by the length of data 'T' as 

Af ^ ^ (uncertainty principle) 


4 


This is similar to a problem encountered in processing speech 
signals where the length of the data stream is limited due to 
the time varying nature of the vocal signals. In fact, it 
parallels the problem clcsely if we consider large cunplitude 
motion where the dynaunics are essentially non-linear and a 
time varying small increment linear model for the airplane 
auid wind is needed. 

The method used there is to model the signal using co- 
variance/correlation mev’-jds digitally. This is based on the 
least squares norm minimization and is referred to as MEM or 
maximum entropy method. 

The method relies on statistical covariance/correlation 
approach to estimating the poles of a digital filter using 
short streams of data and has the desirable property that it 
eliminates the restriction imposed by the uncertainty 
principle by effectively "making more data" using an analytical 
model of the system. Of course, this method has the obvious 
drawback that if the generated model is erroneous, then the 
data generated using it is meaningless. However, the method 
has been successfully used in the speech case and is presently 
being applied to speech recognition. 

Further, the computational and memory requirements are 
quite modest and ein estimate of the accuracy of the model is 
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based on the residual energy. It can be used to decide when 
the model order is high enough if a suitable criteria on the 
size of the residual is available. This was the method chosen 
for use. 


Application to the 
magnetic balance 


The procedure followed was to adjust the model 


suspension to be very stiff in all degrees of freedom except 
pitch. The pitch position response was then stored during a 
run and analyzed using the MINC. For this case all displace- 


ments except pitch are negligible and the equation of motion 


for a second order system with a natural frequency and a 


damping ratio C is 


d^0 . d0 , 2q . 

at * “n ® ■ “ • • • ^ 

dt 

where M is the forcing moment and 0 is the angle of the body 
cixis from the tunnel axis. See Figure 1, 


Z 



Equation 1 in a more transparent form is 


^ % K * V 


dt 


P(t) + M^(t) 


la 


M, 


D„ = 


J 

X 

P(t) 


Aerodynamic moment function 
Magnetic damping 

Moment of inertia about centre of rotation 
Stiffness due to magnetic suspension system 

Pitch input driving function 


The aerodynamic moment can be linearized as follows. 

An Mq term that is due to the lift force due to an angle 
of attack. An M^ term that is due to rotation of any point 
along the body and an M^ term that comes about du«» to 
motion of points along the body (except at centre of rota- 
tion) in the z direction. The total damping M is a com- 
« • 

bination of a and G terms. When the centre of mass is fixed 
as it is in a wind tunnel, G=a, z=0. Hence la becomes 


d^e 


dG 


J ^ ^ ‘V“’ ^ 


In the Laplace domain, neglecting initial conditions. 


0(»){js^ + (Dj^^-M)s + (Kj^ -Mq)}= P(s) 


|i|i = . (K„^-Mg) 
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If P(t) is vhite noise (or an impulse), P (w) is rlat. 
Hence equation 2 can be written more simply in the fre- 
quency plane as 


constant 

0(5) 


= JWp + 




The bandwidth at half power point of the corresponding 
power spectrum can be shown to be (3) 

AWp ^ 2 


where o * frequency in pitch. Using this we have 


AWp = 


"p 


'M. 


M=0 with the wind off and we can measure 


J is gotten 

U 

either from the geometry of the model and its mass or using 

a step input and measuring the natural frequency 'see 4 for 
details) . 


M = J(Aw_ 

^off 


AWp ) 

on 
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M 
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<* 

Hence 

"i “6 


jpu. 
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Now for pure plunging motion in the z direction we 

have 


m 



\3t* 


V 


= L(t) + I^(t) 


.6 


Again using the same argument as above, if L(t) 
is white noise and with o = angle of 

-1 z 

attack - tan (~) . 

Nose up to the wind ^ positive angle of attack 
dictating the sign to be +. 


t 



Figure 2 Origin of lift from a. 


* 

See Etkin. Dynaroxcs of Atmospheric Flight. Wiley, 1972. 
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There is also a lift force due to a. Now for small 


z 


u 


tan'^ (i ) 




- • A 2 

and <» ® ^ 
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Hence, equation 6 becomes 


m 




dt‘ 




and again if L(t) is white, then we can experimentally 
determine 

1 j * L 

(lu - and (D„ - 

"« "l “» 

by measuring the bandwidth at 1/2 power points as 


A(i), 


<°M, - S;> 

I# ** 

(m-^) 
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and » 

n. 


'M. 


(m - 
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or 
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Having gotten and we can relate them to 




1 — 



L< 

a 
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Hence, C, 


.(C 


a 


M: 


+ <=M.> 


are determined. 


In principle and C_ can be determined from a 
n« It* 

a a 

plunging oscillation by measuring the currents and forces 
directly. However, two effects Tiake this more difficult 
tnan using natural resonances. First, because of the large 
inertial loads, forces must be determined from differences 
in large numbers which have both phase and £unplitude. Second, 
a small offset between the center of mass and the center of 
magnetization produces a moment as a result of vertical 
acceleration which must also be subtracted out* Further ex** 
perience is needed to determine if these methods can be used 
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effectively and if these techniques can be used in the side 
slip planes. 

The first step in incorporating Digital control of the 
magnetic balance system was to provide force and position 

I 

inputs to the computer and position set signals to the balance. 
Figure 3 shows the block diagram for this interface. 


COMPUTER INTERFACE 



Figure 3 Computer Balance Interface. 
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A Digital Equipment Corporation minicomputer, the 
NINC Il/:23, was used for the estimation of the stability 
derivatives* The E.P*S. signals from the balance provide 
six positions which, together with the six magnet currents, 
were sampled at 385 Hz simultaneously and recorded on an 
8-inch floppy disk. 

Two Data Translation analog-to-digital boards were 
used — a 16 channel single ended board that was used 'pseudo 
differentially' (the common side of the instrumentation amp 
was referenced to the magnetic baleuice ground) and an eight 
channel differential input board. 

The signal levels being very low (typically < + IV) 
the boards were set up to operate with a full-scale range 
of IV. The 14 bit boards hence had a resolution of 

or approxintately 0.1 mV. 

The cables were shielded, twisted pair and there did 
not seem to be any significant degradation of the signals 
as seen at the computer. 

An extennive library of support programs was written for 

1. Input/output of data 

2. Display of data 

3. Plots of data, and 

4. Processing of data. 

Appendix A lists some of the programs that were felt to 
be useful in the understanding of the processing technique. 
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Most of the actual processing programs are interactive 
and require very little effort on the part of the experimenter* 
A concise procedure for the estimation is as follows: 

I With the model hanging at the desired position, the 
Program *Getdta' - gets in the positions and currents and 
saves them on a user specified file. 

This process is repeat-^ for different wind speeds. At 
each operating point four data files are recorded: 

a. A wind-off noise off data 

b. A wind-on noise off data 

c. A wind-on noise on data 

d. A wind-off noise on data 

'nois?' is the white noise source that is input to the de- 
sired degree of freedom. 

II The Progrcun 'Zerome' - (Zero means) is next run to 
take out the dc characteristics, and the operating point dc 
offsets for each file. 

Ill The Program 'DISPLA' - displays any given data file 
and plots a graph of any given parameter. It also has the 
capability of doing a Fourier transform of the data and 
displaying the FFT. 

This program also plots the graph on a line printer if 
desired. This is basically a prograun to evaluate the data 
and make quick 'by sight' decisions as to the validity of 
the results. 
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IV Program 'FILTER' filters out the dynamics of the 
magnetic balance. 

V The program 'Model' finds the performance function 

I 

of the system using a l^‘ast squares error criteria. This 
is again an interactive program that uses a guess of the 
model order. The residual error energy is monitored and 
the model order chosen that has a reasonably small residual 


energy. 

'Model* also performs a Fourier transform of the 
pe^rformance function which can be plotted. 

^ The half power bandwidth of this resulting spectrum is 
measured and the appropriate parameters calculated. To 
calculate the initial balance parameters and the moment of 
inertia: 


The program 'IGTDTA* - (Impulse and get data) may be 
used. This pulses the model with a small ( 1^) impulse 

and records the resulting damped oscillations. Using 
methods described in (3), the J and the natural frequency 
are evaluated. 



Figure 4 Size of blocks correspond approximately 
to computation time. 
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Wind tunnel tests and results 
including error analysis 

Wind Tunnel Testa - The model was aligned with the tunnel 
axis in yaw and slip and lift and centered in drag* It was 
then set up at a 2-degree pitch angle. 

White noise at 2 intensities was used as an input 
+ 0.125V maximum and + IV maximum which corresponded approxi- 
mately to about 0.05 degrees and 0.2 degrees. Later, repeat 
runs were made at 0.05, 0.1, 0.2 and 0.4 degrees. Data was 
taken for four conditions: 

1. Wind off noise off 

2. Wind off noise on 

3. Wind on noise on, and 

4. Wind on noise off. 

These are shown in Figures 5-8. 

Except in one region, the spectral power did not change 
when the wind or the white noise was turned on, as can be 
seen from the figures. The magnetic suspension feedback 
system poles are distinctly visible with the one varying peak 
corresponding to the model. 

Since the amplitudes of the various other peaks are com- 
parable to the one due to the model, the analytical model 
formed will not be a very accurate measure of the airplane 
model. The data was first filtered using a digital filter 
of length 64 using the Ramez exchange program (2). 


16 










Figure 7 FFT of wind on noise off data. Note the 

reduction of model unsteadiness I (around 11 Hz) 







A pass band was used extending to 18 Hz and a stop band 
extending from 22 Hz. The pass band ripple was designed to 
about 0,65 db and stop band ripple was -20 db. 

The choice of this filter was made arbitrarily. If 

I 

better accuracy in the parameters is desired, the length cam 
be increased at the cost of computation time. The effect of 
the filter was dramatic, as can be seen from Figures 8 and 9. 

It is worth pointing out that the remaining points 
(past 22 Hz) cannot be just chopped off because of the 
ringing that results between samples (Gibbs phenomenon) . 

It is also worth pointing out that the modified signal 
is still only a small signal length (original sample time) 
and hence, still governed by the 'uncertainty principle'. 

The signal was modelled using* an autocovariance method 
and the stationarity of the residual energy v;hen model 
order five was reached talcen to be a. sign of good fit. 

The Fourier transform of the model position is shown in 
Figures 10 and 11 for the wind off and on respectively. 

The excellent resolution in frequency is noteworthy. 
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Saunple Calculatloria 


The model used had the following properties: 



I B 12.5 cms* B 0.12Sm 

p « 1.19 — at 0 altitude, T * 20^C 
meter'* 

m B 345 gms. 

12-3 2 

j pu^ B 4.92 X 10 Newtons/meter 

1 2 

^ PU^SO B 0.06 Newton meters 
1 2 

5 PU^S * 2.41 Newtons 

= 300 ft/sec B 90. 9m/ sec 
O » 0.0254 meters (1 inch) 

Figure 14 is the white noise, as seen on a scope and Figure 15 
is a pitch position response, as seen on the computer. 
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Figure 14 


Typical white noise input. 
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Figure 15 


Typical model response in pitch. 
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Results 


Results of the measurements are summarized below in 
Table 1. Runs A and B were made first. The extensive re- 
pairs were made on the power supplies « the controls were 
readjusted, and a week later Runs C, D, E and F were made. 
From examining Table 1 it is apparent that resolution is 
good but scatter from Run to Run is large, particularly with 
respect to the natural frequency. 
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Table 1 


Preliminary Pitch Damping Data for Ogive Cylinder 
Model at 300 ft/sec, a=*2®, Rejj«0.5xl0^ 


(a) 

Noise ^ 

Amplitude n Am 

Degrees radians 


Damping 

Ratio 


a q 
Dimensionless 


f t/s< 


A 

• 

o 

i/I 

86.8 

20.7 

0.11 

- 

0 

A 

.05 

79.7 

14.2 

0.09 

+ 154.6 

300 

B 

.2 

80.3 

1.8 

0.011 

- 

0 

B 

.2 

73.8 

4.7 

0.032 

- 71.6 

300 


Repeated 1 week later 

after balance 

repair and tuning 

C 

.4 

87.4 

12.4 

0.071 

- 

.1 

0 

C 

.4 

76.2 

18.3 

0.12 

- 143.1 

300 

D 

.2 

93.9 

12.4 

0.066 

- 

- 

0 

,2 

92.1 

t 

20 

0.109 

- 186.1 

300 

E 

.1 

111 

42.5 

0.19 

- 

- 

E 

.1 

86.8 

i 47.2 

0.272 

- 114.5 

300 : 

* 

F 

.05 

119!8 

46.1 

0.192 

- 


F 

.05 

115.1 

40.2 

0.175 

+ 150.3 

300 I 

% 

This 

is a conventional 

*C 

notation with » — and a 

M* j • M _ 

3C 

M 

) * i 

9 (50-) 

It is 

obtained 

from the 

values 

of Reference 4 

by dividing by 

d 

50. • 


This indicates that the control system should have been ad- 
justed to produce the same wind off natural frequency for 
each data point. Also, a larger data sequence should be 
tried to reduce scatter. 

The negative damping (positive coefficient) observed at 
.05^ amplitude is probably a result of system noise. Since 
noise tends to broaden the response curve, it is surmised 
that at very low excitation the decrease in noise resulting 
from damping on the wind offsets the broadening caused by 
aerodynamic damping. At higher excitations where the only 
power points are well above the noise level, the effect is 
not important. This suggests that excitation should be .2 
to .4 degrees. 

An estimate of the relative error can be made as follows: 






,5m 


,5A(i) 




M' 


fiA/P, 


a7d 


/ 


,0.005, 

TTI5- 


,2x0.5, 


i-r s) + 'To(r^ '^75 * 'lyrs' 


,0.05, 
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e = 0.041; i.e., - 4% error. 


These are typical values. 
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This holds for repeated experiments at the same 
eunplitudes. 

The scatter in the three experimental values of runs 
D,E,F (dynamic stability derivatives “143,-186,114) may be re 

la ted to the shift in stiffness in runs. The change in value 
with €implitude may be due to an additional error made in th<? 
fit due to the S/N being different for the three cases but 
further experiments should clarify these issues. 

As with any new technique of magnetic balance experimenta- 
tion, a sizeable base of experience is needed to establish the 
best way to obtain a given type of data. 


Summary and conclusions 

A technique using random noise excitation for estimating 
dynamic stability derivatives using the magnetic suspension 
system has been implemented. The method is not sensitive to 
drift and noise in the system. It also has the advantage of 
being simple and computationally not demanding. A 64K Byte 
MINC was adequate. Typical turnaround times from the start 
of the experiment to the evaluation of dynaunic stability 
derivatives is of the order of an hour. 

Theoretical calculations based on linear theory and 
assuming an = 2.75 inches from the nose, give a value 

of “286 for for this body. 
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Theoretical calculation of C,^ gives a value of 1.477. 

“a 

The experimental values obtained by measuring natural fre> 
quency are 3.32, 0.9 and 5.94 corresponding to 0.1^, 0.2^ and 
0.4° noise. 

Since is computed directly from the change in stiff- 
en 

ness wind on, the scatter is probably due to the change in 
wind off stiffness from run to run.v this can be eliminated 
data scatter should be greatly reduced. 

More work is needed to determine the cause cf this 
scatter and whether static derivative estimation using this 
method can be competitive to the method of measuring currents 
and inferring the forces and moments acting on the body. 

Some interesting facts about the balance can be seen 
from Figures 5 to 9 with the wind on, the balance is 
very quiet. Compare this power (Figure 5 ) to the spectrum 
with the wind on (Figure 7 This says that tunnel un- 
steadiness alone is not enough to excite the frequencies in 
the model for reasonable resolution. 

It is reassuring to note that the system model derived 

2 2 

by the white noise method (s +12.44s+87.6 ) is very close to 
that obtained by Peter Luh using forced oscillation techniques 
(s^+12s+85.21^) (5). 

Two hundred and fifty-six data samples were taken at 
385 Hz, which corresponds to approximately 1 second of data. 
This sampling rate for twelve channels is more than propor- 
tionately lower than the 15 KHz possible for only one channel. 



It may be possible to increase this to nearer 1 KHz with soft- 
ware changes. However, increasing the seunple rate well above 
the natural frequencies of the signal has adverse effects. 

It tends to cluster the poles right on the unit circle and the 
resolution of poles close together is reduced. 385 Hz is more 
than adequate since the highest frequency in the system seems 
to be around 30 Hz. 

The simultaneous acquisition of data from all channels has 
the added advantage of being very useful for the evaluation of 
cross couplings between different degrees of freedom. 

Essentially the seuae procedure may be used to estimate the 
effect of an input in one degree of freedom on another. 

Static stability derivatives may be measured simultaneously 
to the dynamic ones by computation of the mean value of the 
signal with the wind on compared to the mean value with the 
wind off. However, this method has not been compared with 
direct force measurement for accuracy. 

The pareuneter estimation programs can be incorporated as 
a subroutine in another program that analytically computes the 
half power bandwidth from the system performance function and 
hence, directly outputs the dynamic stability derivatives. 
However, a direct involvement in the mechanics of the procedure 
was retained in order to make sure that the nvunbers made sense. 

The 64K Byte memory limitation on the MINC restricted 
taking more than 256 samples at a time.. An additional 64K 
which is on order will eliminate this limitation. Ideally, 
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the number of samples should be approximately 100 times model 
order for a good estimate of the performance function. 

The length 64 FIR filter has a stop band ripple of -20 dB 

and pass band ripple of 0.65 dB. This was arbitrarily chosen 

* 

and the accuracy of the derivatives may be improved by in- 
creasing the length of the filter. This will increase computa- 
tion time, however. 

An important thing to keep in mind is that different model 
dynamics require a change in the filter characteristics (see 
Appendix D) . 

More experiments need to be done with models whose 
characteristics are either known or analytically determined 
and the results evaluated. Once this is done the system is 
ready for use as an easy to use and extremely versatile balance 
that can perform experiments until now not possible. 
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Appendix B 


This is the channel table for the saSnetic wind tunnel coneuter 
interface* 


jdat^r file] /data filel 

w arrau J |u arras J 

cle ^ Sinsle Ended cl»N,yCle E 


ch 

1 

2 

3 

d 


1 

2 

3 

d 

5 


pi^d Position' 
Slie Position 
Yaw Position 
Lift Position 
Pitch Position 


0 

1 

2 

3 

d 

S 


n 

2 

3 
d 

s 

• • 

4 


-4 — 7 


Differential Xnrut 
Lift Current* 

Yaw Current 
Pitch Current 
Brae Current 
Side force Current 
HaSnetizine Current 
Xneut Sienal ~ * 


Data files and their interr retation — * 


UNDFFn.DAT 

UNDFNn.DAT 

WNDNFn.DAT 

UNDNNn.DAT 

ZUDFFn.DAT 

ZUDFNn.DAT 

ZUDNFn.DAT 

ZUDNNn.DAT 

FWDFFn.DAT 
FUDFNn.DAT 
FUDNFn.DAT 
FUMNn.DAT ’ 


Raw data with wind off noise off 
Raw data with wind off noise on 
Raw data with wind on noise off 
Raw data with wind on noise on 

Zero wean wind off noise off data 
Zero a>ean wind off noise on data 
Zero aiean wind on noise off data 
Zero eean wind on noise on date 


Filtered 
Fi Itcred 
Filtered 
Filtered 


wind off noise off data 
wind off noise on data 
wind on noise off data 
wind on noise on data 


(The n denotes the run nua>ber) 


Gauthaa Raaohalli 
H* I * T> * Caiibridie 
8-eiaw-8i 
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C THIS PROGRAM GETS XN THE « CURRENTS* S POSITIONS AND THE INPUT 
CFROH THE MAGNETIC DALANCE AND SAVES THEN ON DISK ON A FILE 
CTHAT YOU ASK IT TO 

C Gauthan Raaohalli 

C H*1«T* tCaabridMc ’ 

C 

PROGRAM GETDTA 
LOGICAL*! LINiDLANK 
DIMENSION LIN<10)*DLANK<10) 

INTEGER ItJfK 

REAL X<2S*iS)»Y<2S6i7)»SUM*MEAN 

INTEGER XX(2S*I 

N-254 

DATA X*Y/1280*0«»1792*0«/ 

• WR1TEC5»S00) 

500 FORMAT (* ENTER 10 BLANKS^ 

RCAD(S*S01« BLANK 

501 FORMAT (lOAl) * “ 

1 DO SSS 1-1*10 

tIN(I)-BLANK(I) I 

SSS • CONTINUE 

CALL A2liM(0*0tS*N»X*l.) 

CALL A2DM(l*0i7»N*Y»l.) * 

CALL SCAL(Y(1*7)*1X«N) 

SUM-0 ^ 

DO 200 I-1*N 
SUH-SUM^Y(1>7) 

200 CONTINUE 

MEAN«SUM/FLOAT<N> > 

CALL CLEAR 
CALL GRAPH(N*1X> 

URl TECS* 990) MEAN 

999 FORMAT<' Naan Valua-* *F8*S* ' Tha data diaalsuad i« iha whita noiaa 
X in^ui*) 

2 URITECS*998) 

998 FORMAT<' Do uou uant this data wriitan to a filaT u-l*Raad ao 
Xrf-O*ool t--2' ) 

READ(S*997) KR 
lF(KR.EQ.-2> CO TO 1000 
IFiKR.NE.l) GO TO 1 
CALL CLEAR 
UR1TE<S*99A) ' 

997 FORMAT < 12) 

99A FORMAT!'’ Entar tha naaa of tha file to ba saved as') 

READ(5*SS7) LIN 
SS7 FORMAT! lOAl ) 

OPE N ( UN 1 T- 1 • NAME-L 1 N* FORM- ' UNFORMATTED' • TYPE- ' NEW ' ) 

DO 1200 I-1*N 
WRITE!!) !X(!t J)»J-1*5) > 

URI1E!1) !Y!1* J)*J«1*7) 

1200 CONTINUE 

CLOSE !UN1T-1) 

CO TO 2 

1000 CALL CLEAR 
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C 

C This preiras takas an input fila and filtars a divan channal 
C and writas ^ha fila bach as a filtarad \*ila* Tha prodraa usas 
C a 44 lanath finita ispulsa rasponsa filtar with a pass band 
C that awtands to 18Hz and stop band that starts at 2?H2« This 
C affpctivalw filtars out tha nadnatic balanea swstao faadbacK 
C affacts froa tha sidnal and what raaains is • pura sidnal 
C dua to tha airplana dwnaaicsl 
C ... 

C Cauthas Raaohalli 

C 8-aav-81 

C H* X«T* tCaabridda 

C 

Prodraa filtar 

diaansion LIN(10>»L0UT<10)»BLANK(10> 

LOOlCALdl LINtLOUTfBLANKf rtb . 

Raal K<2S4fS>rs(2S4»7>»h(44)trdua<7) 
coaploK mkC 1024> »hh<10241 
Intadar i » Jtkf ichnltflad 
I data b/lHP/ 

’ OPEH<UNlT-l»NAMC«'FILTER.BAT'tTYPE-'OLO*r' * 

READ<l»$Mh(i}*i•l«32^ 

I closa ( uni t*l > 
do 8081 i-l»32 
h(it32)>hC33-i) 

8081 continua 
do 8082 i«li44 

hh<i)>h(i) 
writa(4f8083> i»h<i) 

8083 foraati* hC^fI2»M « '»C20.f) 

8082 continua 

do 8084 i-4S*1024 
hh(i)>0* 

8084 continua 

call f ouraaihht 1024»~1 I 
URITE(S»?S0) 

850 FORMAT (' ENTER 10 BLANKS PLEASED 

READ(StSOl) BLANK 
501 FORMAT <10A1) 

1 DO 1000 i«l»10 

L1N<1 )>»LANK(i> 

LOUTf i)>BLANK(i) 

1000 continua 

writa(S»99?> 

988 foraatC* antar tha naaa of tha input fila*) 
raad(S*888) LIN 
888 foraaKlOAl) 

wrltc(5t8S8> LIN 
writa(Sf27U) * 

2714 foraati* antar tha channel to ba filtcralKcd r5) *> 
raad(S.2717)r>lchnl 

2717 foraat<AltI2) 
writaCSf 2718)r»ichnl 

2718 foraati* filter channal *«Altl2**T*) 
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♦ . . i ^ * * 

r«ad(8»271f) kr 
2710 for»at(I2> 

if (kr.aacS) do to 1 
OSB foroati* fllo riaot««f lOAl) 

1040 continuo 
flao-0 

if(r.na.b» ftno*! 

0PCN< UNlT-ltNAHC-LINtFORM-^UNFORHATTED'iTYPE-^OLG'l 

do 1100 i«l»2S4 " ‘ 

raad(l) (K<it J) t J*1 »S> 
roadUl (o<it J) » J«lt7> 

1100 continuo 

clooo(unit*l> 

if (f lad«to«l)do to 1200 

do IISO i-lt2S4 

kmC i)«M(itichnl) 

IISO continuo 

00 TO 1250 

1200 do 1250 i-lf254 

»K(i><*o(liichnl) I 

1250 continuo ' 

do 1270 i-2S7tl024 

KK<i)*0*0 t 

1270 continue 

call fouros<x)(tl024t**l> 
do 1280 i«ltl024 
MM(i)-KKCi}thhCi) 

1280 continuo 

call f ouroaixxf 1024* 1 1 
if<flao.ao«n do to 1310 
do 1290 i«lf2S4 
Kfi»iehM}-Roal(xM(i>) 

1290 continuo 

do to 1320 

1310 do 1320 i«l»2S4 

v< i f ichnl )«Roal(xx( i 1 1 
1320 continuo 

997 foraiat(/t' ontor tho nano of tho outaut filo'l 

957 FORHATi' iAaut f ilo«^ tlOAlf * ouiaui«^«10All 

«ritt(S*997> 
road(5t998) LOUT 

urito(S»957) LlNtLOUT . . 

OPEN ( UN I T «2 1 NAHE-LOUT t FORH- ^ UNFORHAT TEIl * t T YPE> 'NEW' > 
do 1500 i-l»2S4 

writa(2) Cx<itJl*J-l»S> 

Mrito(2) <«<if J) t J*1«7I 
1500 continuo 

clofta(unit«2) 

«iritoi5»990) 

990 foroat(/t' Do oou havo anu aioro filoo to work on? o-l'l 

raVd(5«989) kr . . 

989 foroat(i2> 

if<kr*oo*ll do to 1 
fttOO 
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C 

c 

C and 
C Th« 
C and 
C for 
C 
C 

c 

c 


Thi« rreirao raada In a inrui filo < nadnoife balanca data! 
rarforaa a covarlanca ntihod of aidnal aodalling 
rrodraa ia iniaracilva and ask* for tha fila nano 
iha landih of data and aodol ordar* li loora back 
a Chanda in nodal ordar* ‘ ' 

Oauthan ftanohalli 

3 >aaa-ai 

H*1 *T* fCanbdridda 

Frodran Nodal 
' lodicaldl LIN»blank»r »b 
diaanaion LIN(10)tblank<10> 

doubla rraeiaion X<2S4)tA<2\)*K!:C21)fdd<300}»A2C2ll 
doubla rraeiaion ALPHA ' * 

r«al rdua(7}*w(S12)»KMCS12)>PX 
iniadar nntotar t i i Jtkt 1X(512) tff * iehnl 

conrlaM z(1024>fH 

bATA b;d/lHPtS12kO*/ 
ff«12 

PI-ATAK(1.)«4. * . . 

U«CHPLX<0.fPn 
arilalStlS) b 
f oraaK IxrAl ) 
no 17 i-lt25d 
X(i)>0. 
coniinua 
do It 1-1*21 
RC(i)-0« 

A(i)-0. 

coniinua 

ALPHA-0. 

nn-29A 

writa(S»SOO> 

foraaii* Cntar a blank lina rlaaaa') 

raad<5tS01> blank 

foraaUlOAl) 

foraaii' Cntar iha naaa of tha fila to ba raad'l 

do S i-l>10 

LlN(i>-blank(i) 

coniinua 

uritaiSfttt) 

raadiStttO) LIN 

foraaUlOAl) 

uritatS*tt7) LIN 

foraati'l Fila naaa la *»10A1) 

writaC9itt4) , 

foraatt'lTo diarlau a rarticular raraaatar tura tha arrr no^l 
writaiSfttS) 

foraatl/** Pl»drai roaition* »T40* 'Cl-1 ift currant*} 
wrna(S«tt4> 

feraatC* P2-alir roaition* «T40* *C2-uau currant*) 
writa(S*ff3) 

foraatC* PSa^uau roaition* * T40* *C3«rltch currant*) 
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* urit«(5t992> 

' 992 forkaK* P4«lift position' tT40t 'C4«dratf current')* 

uriie(5»991> . , _ - - 

‘ 991 foraaK' P5*pitch position' tT40t 'CS«slip current') 
write(S»990> 

99y foraat(T40f 'C6«asSnetizinn current')" * " ‘ 

write(5»989) 

989 foraat(T40t'C7«input sidnal') 

< read<5»988) rtlchnl 

r 988 foraat(Al»12) 

v'iie<Sf967) rtichnl 

987 .9 aatC'l't' display 'tAlrI2t 'data?* (-*1 to read file adain)') 

r «dtS»98&) kr 
it(kr<^ e*-l) do to i 

‘ •"call clear - - . 

if(r«ee«b) do to 1190 

0PCN(UNIT«ltNAHE«LXNtF0RH«'UNF0RHATTE8' vTYPE>'0LD') 

■■ do 1100 i«ltnn 

read(l) <rdua( J) » J«lt5) 

read(l) <rdua( J) * J«lt7) • ( 

y(i)«rduB(ichnl) , . i . 

X< i)^rdua(ichnl) 

1100 continue ) 

close(unit*l> 
do to 1200 

1190 OPEN ( UN I T « 1 » NAHE-L I N t F0RH« ' UNFORMATTED ' t T YPE- ' OLD ' ) 

do 1250 i-l»nn 

readCl) (rduaC J) » J«lrS) 
y(i>«rdua(ichnl} 

X<i)*rdua(ichnl> 
readCl) <rdua< J> t J*lt7) 

1250 continue 

close(unit*l> 

1200 continue 
1280 call scal(ytXX>2S0) 

call draph(250tIX) 

1275 urite(5»985) 

985 foraati' ret to read neu dtal2 aiodel sidnalli**! to ouit)') 
read<Sr98&) kr 

* 984 foraat<I2) 

call clear _ 

if(kr.ea*'l) do to 3000 ' 
if(kr»ea.2> do to 2500 
do to 1 

, 2500 write<5t969) 

. 949 foraati' Make sure the line printer is set up ridhtl') 
urite(S»948) 

948 foraati' The vertical lpi«4»foralendth«44') 
read(5»947) kr 
947 foraat(I2> 

SS7 write(4t777> ff 

Hrite(4»970) LINtrtichnl 

970 foraat(///t' Dataf ile-' tlOAl >T30t 'Pos/Cur-' t A1 tTSOt 'Chnl»' 1 12 
777 foraaUAl) 
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979 

urite(6*979> 


FuKniiTfrf' . 





'*write(5>950) 


950 

forsati* oodel 

order*?*) 


read(5»951) ■ 

• .... 

951 

forsat<I2) 


952 

write(6t952> s 
f ornate* aiodel 

order** tl3) 

953 

write(5t953> 

fornate* lensth of siSnal>T*> 


' read(5*954) n 

* * — 

954 

fornateit) 



call covar(n*Xf 

•»A*ALPHA»RC) 


*- • 

«ritc<At9SS) n 

955 fomaK/t' length of th* si8nal«*tX4) 

~ ~uritef6»9S4) " ' 

956 for»ai(/»T22» 'hodel coefficients^ »TS0t'rc.'f lection coefficient 
do 350 i«lfoe 

wriie(6i9S7) (itACi>»RCCi)ti-lt>e)‘ 

957 for»at(/fI3f 17XtE20*10»10XfE20.10) 

350 continue 

write(6i9S8) ALPHA 

958 foroai</t' Residual enerd«<K' t E20* 10) 
if<».)t.7> 90 to 8080 

- wriie(4»777) ff 

8080 write(6*979) 

iduo«2»i.4256 
n>idus 

do 937 i-lto / 

V9<i>-0. 

' ««( iduo-iiDaO* 

937 continue 

do 959 i-lt254 
MV<i4»)«X(i) 

959 continue 

call auioinr yyrBi*A2t ALPHA* RC> 
urite<6'?955) n 
wriie(A*956) 
do 800 

writc(4*957> ( i * A2C1 > *RCf i > t i-1 »t.e) 

800 continue 

wrile(A»958> ALPHA 

URITE(6»979> 

write<5i555) 

555 forsstC' ret to chande oodel order! 1 to read new filel 
X 2 to set ffts of oodelst-1 to ouitlM 
read(S*98A) kr 
if(kr*ne«2) So to 8092 ; 

write(5t7070) . 

7070 fornati* which nodel do wou want ffts for covar*0*auto«l* > 
read(St98A) kr 
flas«0« 
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ifCkr*ea*l> fladsl* 
tiriie(5t8089) 

8089 fomatC'^ co»^uLind ffts') 

do to 8082 

* G do 8082 

C z(i>«A(i> 

8082 continuo 

C do 7072 

C z(i)-A2Ci> 

7072 continuo 

C do 8083 i*»^ti»1024 

C z(i)*0« 

8083 continuo 

C coll foureofz* 1024? ~1) 
do 707 i*>l?S12 
z(i)«i.O 
do 702 J«2fOO 

z(i)«z(i)4A(J>«CEXP<-<U«FL0AT<i-U/2048*)«FL0AT(J-ll) 

702 continuo ■ 

z<i)>l*/z(i> 1 

707 continuo 

do 8084 i-lfS12 > 

i >»cabs(zC i) ) 

8084 continuo 
call clear 

call 9cal(o?IX?S12) 
call drarh<S12tIX) 
call shadow 

8091 write(Sf8085> 

8085 foroati' rot to loorl ^ to rrint riot I '-‘1 to ouit') 
road<S?980) kr 

if(kr»no«2) do to 8092 | 
do 8080 i«l?S12 

xx(i)-(floatii-l)/4090.)«38S. 

8080 continuo • 

writo<S?8087) 

8087 f ornate ' entor no of roLnts to bo plotted') 
read<S»8088> nrts 

8088 foraat(13) 
writo(0?777) ft 
write(0?970> lin?rtichnl 
urito(0?979) 
writo<0?9S2> • 
urite(0f979) 
writo<S?8090) 

*8090 foroate' oako sure printer is sot for 12 1p1 and 132 Ipf') 
pause 

call plot(MX?wtnpts) 
call clear * 
so to 8091 

8092 call clear 

if (kr)3000fSS7?127S 
3000 call clear 

stop 44 
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LIST OF SYMBOLS . 

% 

Magnetic damping in lift 

Magnetic stiffness in lift 

Magnetic damping in pitch 

Magnetic stiffness in pitch 

Damping ratio 
Frequency in radians/sec 

Moment of inertia about the center of gravity 
Free stream velocity in the -x direction 

Pitch angle with respect to the x axis 

Moment acting on the body 

Angle of attack 

Dehsity of air 

Mass of body 

Characteristic area 

Characteristic length 

By definition is 





t 


L = aerodynamic lift force on the body. 




C, =* 


C, * 


3 Cl 


do 

3 C, 


a& 


'M 


I pu.*SD 



3e 



!_5i 


note, this DEFENITIOH LEADS TO UNIT§ ,0? SECS. FOR 


V 
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Appendix D 


FILTER DESIGN GUIDE 

The Remez exchange algorithm needs the following 
inputs: 

Card 1 N, ITYPE, NBANDS, MGRID 

N - order of filter* (This is an inter- 
active procedure where you keep increasing N \mtil the 
desired attenuation in the stop band and pass band ripple 
is achieved.) 

N»l: Low pass high pass or band pass 

NBANDS =: Number of bands, eg., a low pass 
filter has two bands. 

NGRIO: Grid density always set to sixteen. 

For example, a length fifteen low pass filter would have 
the following first card: 15, 1, 2« 16 

Card 2 FILOH, FlHIGH, F2L0W, F2HIGH 

FNBANDS LOW, FNBANDS HIGH 

This card -specifies the edges of the bands. There are 
MBANOS as specified on Card 1. These band edges are entered 
as real numbers (not more than four per line) and normalized 
to (0.0, 0.05) interval corresponding to the interval (0-1/2T) 
where T is the sampling time. In this case 1/2T “ 187 Hz. 
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Thus a low pass filter with pass band to 18 Hz* stop 
band from 22 Hz would have Card 2. 0.0, 0.0467« 0.057# 

O.S 

Card 3 AMPl,.. AMPNBANDS 

This card specified the desired magnitude response in 
each band. For example# a low pass filter would have 

Card 3 1.0# 0.0 

Card 4 HI# ... HNBANDS 

Specifies the ripple weighting in each band. These are 

real numbers. The ripple' ratio is given by 4^/4^# where 

i 

4^ is the maximum derivation allowed in the first band 

and 4^ is the derivation allowed in the ith band. 

Thus for a low pass filter with 4^^ * 0.1# 42 * 0.2# 

this card would be 1.0# 0.5 

The I'esult of running this program is a data file 
called FILTBR.DAT that gives the impulse response. This 
file should then be transferred to the disk where the 
FILTER. SAV program resides (i.e. , the same disk that you 
have your magnetic balance data). The FILTER. SAV program 
reads in the filter characteristics before it operates on 
the magnetic balance data file* 

Detailed instructions for the use of this procedure 
will be found in a user manual that is being planned and 
should be available at the balance. 
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t*t****t**tttt*t***tt*tt*t*t***tt.tttttt*ttt*t*tt*ttt 

* « 

Gautha» Raaohallir 
U91/106f (3-8334> 

H*X«T*f‘ CanbriddSf Ns* 02139 


-••INPUT FILE 


34 1 2 

0.000000000 
1.000000000 
2.000000000 


13 

0.043799999 

0.000000000 

1.000000000 


0.057100002 


0.500000000 


CALLING REHEZ AND RUNNING PROGRAH* DONT INTERRUPT 


DEVIATION 

DEVIATION 

DEVIATION 

DEVIATION 

DEVIATION 

DEVIATION 

DEVIATION 

DEVIATION 


-0.0002315GS 

-0.039413721 

0.132437500 

0.150773944 

-0.155212839 

-0.153539904 

-0*153335235 

-0.153335449 




FINITE IMPULSE RESPONSE <FIR) 
LINEAR PHASE DIGITAL FILTER DESIGN 
REMEZ EXCHANGE ALGORITHM 


BANDPASS FILTER 


FILTER LENGTH « 34 


ttttt IMPULSE RESPONSE ttttt 


H( i> 

m 

-0.73954342E-01 

a 

H< 

34) 

H( 2) 

m 

0.33773184E-02 

8 

Hf 

33) 

H( 3> 

m 

0.4S733V05C-02 

m 

Hf 

32) 

H< 4> 

m 

0.33917830E-02 

X 

Hf 

31) 

Hf S> 

m 

0.8SS73177E-02 

m 

Hf 

30) 

H( 3) 

m 

0.10734314E-01 

m 

Hf 

59) 

M< 7> 

m 

0.1233S142E-01 

m 

Hf 

58) 

H( 8) 

m 

0.13839498E-01 

m 

Hf 

57) 

H( 9) 

m 

0.14070034E-01 

9 

Hf 

53) 

H(10> 

m 

0.13133059E-01 

m 

Hf 

55) 

H<11> 

m 

0.10920897E-01 

9 

Hf 

54) 

H(12> 

m 

0.74827075E-02 

9 

Hf 

53) 

H(13> 

9 

0.29894039E-02 

9 

Hf 

52) 

H<14) 

m 

-0.21979585E-02 

9 

Hf 

51) 

H<1S) 

m 

-0.73348444E-02 

m 

Hf 

50) 

H(13> 

m 

-0.128SS738E-01 

m 

Hf 

49) 

H(17> 

M 

-0.172S37S9E-01 

m 

Hf 

48) 

H(ie> 

m 

-0.19903990E-01 

m 

Hf 

47) 

Hf 19> 

m 

-0.20732537E-01 

9 

HI 

43) 

H<20) 

m 

-0.19143100E-01 

9 

Hf 

45) 

Hf21> 

m 

-0.1484889SE-01 

m 

Hf 

44) 

H<22> 

m 

-0.77870712E-02 

9 

Hf 

43) 
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HC23> 

m 

0.iy024462E-02 

m 

H( 

42> 

H(24> 

m 

0.13866879E-01 

m 

H( 

41) 

H<25) 

m 

0.275661S7E-01 

m 

H( 

40) 

H(26) 

m 

0.42291783E-01 

m 

H( 

39) 

H<27) 

m 

0.S7208348E-C1 

m 

H( 

38) 

H(28> 

m 

0.71434379E-01 

m 

HI 

37) 

H(29) 

m 

0.84089853E-01 

m 

HC 

36) 

H<30> 

m 

0.94392702E-01 

m 

H( 

3S) 

H(31) 

m 

0.10165238ET00 

m 

H( 

34) 

H(32) 

m 

0.1053S892EF00 

m 

H( 

33) 


LOWER BAND EDGE 
UPPER BAND EDGE 
DESIRED VALUE 
WEIGHTING 
DEVIATION 
DEVIATION IN DB 


BAND 1 
0.0000000 
0.0468000 
1.0000000 
2.0000000 
0.0783177 
0.6S493S2 


BAND 2 
0.OS710O0 
O.SOOOOOO 
0.0000000 
1.0000000 
0.1S663S4 
-16.1021996 


EXTREMAL FREQUENCIES— MAXIMA OF THE ERROR CURVE 

0.0117187 0.0273437 0.0400391 0.0i468000 O.OS71000 
0.0649125 0.079560l9 0.0951859 * 0.1108109 0.1264360 
0.1430375 0.1586625 0.1742875 0.1899125 0.2065141 
0.2221391 0. 2377641 0.2543657 0.2699907 0.2856157 
0.3012407 0.3178422 0.3334672 0*3490922 0*3647172 
0.3813188 0.3969438 0*4125688 0.4281938 0.44479S3 
0*4604203 0.4760453 0.4916703 


t*t**tt**ttt*tt*t****t*ttt**t*tttttt*ttttt$t*tt*tttttt**t****tt 

1 

.0 0 0 0 


{ 


. 'SO 



